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Abstract 

A new method to construct event-generators based on next-to-leading order QCD matrix- 
elements and leading-logarithmic parton showers is proposed. Matrix elements of loop di- 
agram as well as those of a tree level can be generated using an automatic system. A 
soft/collinear singularity is treated using a leading- log subtraction method. Higher order 
re-summation of the soft/collinear correction by the parton shower method is combined 
with the NLO matrix-element without any double-counting in this method. 

An example of the event generator for Drell-Yan process is given for demonstrating a 
validity of this method. 
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1 Introduction 



The standard model(SM) of electromagnetic and strong interactions has been well established 
through precision measurements of a large variety of observables in high-energy experiments pQ. 
However, there is still one missing part of the SM to be confirmed by experiments, i.e. the Higgs 
boson. In order to search for this missing part and also search for direct signals beyond the 
SM, the LHC experiments [2j are under construction at CERN. Though LHC has employed a 
proton-proton colliding machine due to their high discovery potential for new (heavy) particles, 
a capability for precision measurements must also be seriously considered. For extracting 
physically meaningful information from a large amount of experimental data, suffered from huge 
QCD backgrounds, a deep understanding of the behavior of QCD backgrounds is essentially 
important. 

The prediction power of a lowest calculation of QCD processes is very limited due to a large 
coupling constant of the strong interaction and an ambiguity of the renormalization energy 
scale. Predictions based on the NLO (or more) matrix elements combined with some all-order 
summation are desirable. Much effort to calculate NNLO corrections j3] is being made for LHC 
experiments. Moreover, in order to suit for precision measurements, a simple correction factor 
(so-called if-factor) for the lowest calculation is clearly not sufficient. One must construct event 
generators including NLO matrix elements with higher order re-summation even for processes 
with multi-particle final states. When one tries to do this task, one may encounter the following 
difficulties: 

• the number of diagram contributing to processes with many particles in a final state is 
very large, 

• many parton-level processes must be combined to construct a proton-(anti-)proton colli- 
sion, 

• numerical instability in the treatment of a collinear singularity may appear, 

• a careful treatment of matching between matrix elements and the re-summation part 
must be required to avoid double counting, and 

• negative weight of events in some part of phase space is not avoidable. 

In order to overcome these difficulties, lots of works have already been donejHl Ej. In this 
report, we propose a new method to construct event-generators with next-to- leading order 
QCD matrix-elements and leading-logarithmic parton showers. Our solution for the above 
problems will be given in this report. 

2 Matrix elements 

A calculation of the matrix elements of hard-scattering processes is not a simple task when 
multiple partons are produced. At the tree level, the GRACE system^, an automatic system 
for generating Feynman diagram and a FORTRAN source-code to evaluate the amplitude, is 
established. The system can in principle treat any number of external particles, and has 
been used for up to six fermions (> in the final state within a practical CPU time. In the 
GRACE system, matrix elements are calculated numerically using a CHANEL 7 j library based on a 
helicity amplitude. CHANEL contains routines to evaluate such things as: wave-functions/spinors 
at external states, interaction vertices, and particle propagator. 

Based on this method, we have already published an event generator of four 6-quark gen- 
eration processes at proton-(anti-)proton collisions at tree-level named GR@PPA_4b[Hj, which 
include all possible tree-level diagram to create four 6-quarks in final state with an interface to 
PYTHIA0. 



2 



For the loop diagram, effective vertices of two- and three-point functions have been prepared|10| 
and implemented in the system. At first amplitudes of the loop diagram are evaluated in d- 
dimensional space-time. The dimensionless coupling constant is introduced as g — > ji d / 2 ~ 2 g. 
After the loop-momentum integration, UV divergences are eliminated by subtracting 1/ejjy + 
je — ln47r where ejjv = 2 — d/2 and 7^ is the Euler's constant. This is so-called MS renor- 
malization scheme [TT], Then we are left with the infra-red divergence which is also regularized 
using the dimensional regularization by setting d = 4 + 2ejr, where eir > 0. 

Expressions of effective vertices after the MS renormalization are summarized in Ap- 
pendix A. For more than three-point diagram, an extended CHANEL is prepared. It can evaluate 
a fermion current including any number of gamma matrices. An amplitude of loop diagram will 
be calculated based on the fermion current connected by boson propagators. When a fermion 
loop is included in the diagram, numerical calculation will be performed after taking its trace. 
A bosonic loop can also be treated by the CHANEL. 



3 Loop integral 

As two- and three-point functions are embedded in the system as effective vertices, loop integrals 
which must be done in each process start from four-point functions. The loop-momentum 
integration of the four-point diagram can be expressed as 

T 1 \ ( 2\ £uv f d d k 1 / -1 \ 

h(pi,P2,Ps,M) = H J {2lT y lAlA2A ^ (!) 

A 1 = k 2 + iO, 

A 2 = (k+ Pl ) 2 + iO, 

A3 = (k + pi +p 2 ) 2 + i0, 

A 4 = (k + Pl + p 2 + p 3 ) 2 + *0, 

Pi + P2 + P3 + Pi = 0, 

where k is a loop momentum and pi's are external momentum taken to be incoming. After a 
momentum integration using an usual procedure of Feynman-parameter method, the integral 
can be represented as 



h(pi,P2,P3,P4) = — jjzxi — Jds,t,ml,ml,rr4,ml), (2) 



r (2 - S ir) ~ t ^ _ 2 _ 2 _ 2 _ 2 
(47T) 

where m| = p| , s = (pi + P2) 2 , and t = (p% + P4) 2 . In general we need tensor integrals 
1/4(5, t; n x , Tiy, n z ) = 

(Attu 2 ^) / dx dy dz n , (3) 

V P / Jo M ( — xzs - y(l - x - y - z)t - iO) 2 ' £lR K ' 

here we set all external particles on-shell (massless). A method to perform this integration 
is well established. The scaler integrations (when = 0) without an infrared 

divergence are given in ^2] an d those with IR divergence are given in |13j . Tensor integrals 
can be obtained from a linear combination of scaler integrals |14| . Even though the method 
is well established, we decided to develop our own method for the tensor integration of an 
IR divergent case, in order to avoid a possible numerical instability due to some cancellation 
among IR divergent scaler integrals. 

The tensor integration can be done analytically and be represented by a finite number of 
terms with Beta and Hypergeometric functions as; 

T , s. _ 1 B(n x + £i R , n y + n z + e IR ) 

J4(s, t; n x , ny, n z ) — — 

s t l - sir 
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AlTfJ, 2 



n , 



■B(l + n z ,n x + n y + e IR ) 



x 2 Fi 1 + n x , n x + n z + eir, 1 + n x + ra„ + n z + 



+ h — 2 L T J n : -B(l + n y ,l + n z + e IR ) 

x 2 Fi^l + l,l + n z + e IR ,l + l + ny + n z + e IR ,-^j , (4) 

s = s + iO, 
t = t + iO, 
i = t-iO, 

u = u + iO = (pi +P3) 2 + iO. 
When the numerator is unity (n x = n y = n z = 0), the scalar integration can be obtained as 

1 B(ei R , ei R ) 



J 4 (s,t; 0,0,0) 



(— 



sir 



s te IR (l - em) 

2F1 I, Sir, 1 + Sir, 



II 



sir 



2F1 1, EiR, 1 + Sir, 



(5) 



A explicit structure of the infra-red singularity can be obtained with Laurent expansion of 
above formula with respect to sir. A full set of formulae and a computer program for the 
numerical calculation can be found in ^f)j 

For more than four-point functions, a reduction formula [TK] is used in the system. 



4 Real emission 

4.1 Soft/Collinear treatment 

At first, the conventional method of phase-space slicing[T7] is used to treat soft and collinear 
singularities. It is explained for the case of an initial-state radiation in this report. The final- 
state radiation can be treated in a similar way. Let's consider two colored partons, whose 
momenta are p^ and p^, scattered into the A-body final state of colorless particles as a Born 
process. For the radiative correction of this process, we must treat processes with one additional 
colored parton, whose momentum is k^, emitted in addition to the Born process. The matrix 




Figure 1: Dalitz plot of (N+l)-body phase space. Qi is defined in Eq.© and is a threshold value to 
separate the visible- and collinear-region. 
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elements of the real emission processes must be integrated in (iV+l)-body phase-space, 
in a (i-dimensional space-time. The space-time dimension is set to be d = 4+2ei R in this report. 
A collinear region in the (N + l)-body phase-space of the final particles in the (i-dimensional 
space-time is defined as 







coll 



a {k»\c% = -(jh - kf < Ql) c n full = *™ +v 



(d) 



(6) 



as shown in figure Q where Q 2 is some cut-off value. It must be sufficiently small so as not to 
be observed experimentally. Final results must be independent of this value. In the collinear 
region, matrix elements can be approximated as 



M 



(d) 

N+l 



M 



N 

£■ 

i=l 



167r . f ^ip^-L 



1 



(7) 



where P(x) is a splitting function for a given parton splitting at a leading-logarithmic order, 
kx is the transverse momentum of the radiated parton, f c is a color factor of a given branching 
and \x the energy scale of the splitting. The CM-energy of the Born process is s = q 2 = sx, 
where s = (pi + P2) 2 ■ Here, we can set d to be four in the matrix element of the Born process, 
since it has no IR singularity. In the same approximation, the phase space is expressed as 



(d) 

N+l 



JV 

i=i 



( k 2 T \ 



SIR ]_ 



16vr 2 r(l + e IR ) \A7rx 2 J 1 



-dxdk 



r- 



(8) 



The total cross section of real-radiation processes in the collinear region can be obtained as 



(2 P ° 1 )(2 P ° 2 )v rel Jn co 



d<£ 



(d) 
N+l 



M 



(d) 
N+l 
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d<S> 



(d) 
N+l 



full 



^"uis 



(d) 

N+l 



M 



{d) 

N+l 



where Q V i S means a visible region of the phase space, such as 







full 
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(9) 
(10) 

(11) 



Then, the collinear cross section can be obtained after the integration with respect to k\ as 



a coll 

O full 



G full @vis > 
1 
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s 



d<S> 
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(12) 
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(14) 

(15) 

(16) 

(17) 
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is the Born cross-section at a CM energy of q 2 = xs in four-dimensional space-time, ji is a 
factorization energy scale, and is a step function. For example in case of a splitting for a 
quark to quark+gluon, a co u after the Sir expansion becomes 



*o(*)g/c 



+ 



-IR 



2L-3 
sir 



71 



+ 



2 / dxa (xs)(f>(x,E IR ) 
Jo 

dxao(xs) 
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(1-x). 
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where 



, \ 1 f as v( \ 1 , a s 1 + x 2 

{x,e IR ) = — fc-z-Hx) = — far- 7: — 3-, 
27T e IR 2tt (1 - x) + 

L = ln(s//i 2 ). 



(18) 

(19) 
(20) 



The IR-divergent terms in the collinear cross section can be canceled out after summing up 
terms from virtual corrections. After combining with virtual corrections, there still exists the 
collinear-divergent term, <p(x,ei R ). This term is thrown away by hand, because it is counted 
in the PDF or PS on the initial partons. Since the MS renormalization scheme is used in this 
calculation, the same scheme must be used in the PDF or PS. 



4.2 Visible jet cross-section 

The last part of the NLO cross-section calculation is a real emission of the additional parton 
into the visible region with exact matrix elements. Those matrix elements of (iV+l)-parton 
production at the final state in four-dimensional space-time can be automatically generated 
using the GRACE system. Cross sections are obtained by integrating those matrix elements 
under the four-dimensional phase-space as 



1 



O 'exact 



(2p° 1 )(2 P ° 2 )v rel Jn v 



(4) 
N+l 



M 



(4) 
N+l 



(21) 



The space-time dimension is set to be four {sir = 0) hereafter, since there is no IR-divergence 
in Q v i s . When the threshold energy (Q 2 ) is set to be sufficiently small, this cross section can be 
larger than the Born cross section due to the large magnitude of the coupling constant (a s ). A 
parton-level calculation has no problem, except for the large higher-order correction. However, 
if one tries to calculate the cross sections of a proton (anti-)proton collision, one has to combine 
the matrix elements with the PDF or PS to construct a proton from partons. The PDF and 
PS include leading-log(LL) terms of the initial-state parton emission. If the matrix element 
is combined with the PDF or PS very naively, a double-counting of these LL-terms must be 
unavoidable. Our proposal to solve this problem is to subtract the LL-terms from the exact 
matrix elements as 



O'LLsub 
fLL(x,s) 



1 



(2p»)(2p°K ei Jn m 
U 2tt P[X) k 2 



(4) 
N+l 



M 



(4) 
N+l 



M^\sx) 



fLL(x,s) 



(22) 
(23) 
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The second term of the integrand is the LL-approximation of the real-emission matrix-elements 
under the collinear condition. There is nothing but 'cr^s' in the last subsection. Then the real- 
emission cross section (a rea i) can be expressed as 



@col ~l~ O exact i 

O full &vis O exact j 

a 'full + O'LLsub- 



&real — &col <? exact > (^4) 

(25) 
(26) 

One can expect a numerically stable calculation using Eq. (j26j) rather than Eq . (|24|) . because 
cancellation between a v i s and a exa ct occurred before the phase-space integration. When the 
threshold value of Q 2 is set to be sufficiently small, the integrand in cr^Lsub is very close to zero, 
because the LL-approximation is very precise around the collinear region. Then, the result of 
integration in (TLLsub is independent of the threshold value of Q 2 . 

This subtraction may distort the experimental observables of additional jet distributions. 
However, a subtracted LL-part will be recovered after adding the PS applied in the Born 
process. The exact matrix element of real-radiation processes gives only non-logarithmic terms 
to the visible distributions in the LL-subtraction method. The relation between the PS and 
the LL-subtraction method is as follows. The LL-approximation cross section in a^Lsub is 



1 



(4) 



(2p?)(2p> re , J Uvis N+1 
dxaQ{xs)V <yl \x) , 



M%\sx) 



fLl(x,s) 



o 2 

ma -j 



dkl 



Ql/il-x) fef 2ir 



(27) 
(28) 
(29) 



where Q max is the maximum value of the transverse momentum. T>^ (x) is just the first term 
of the all-order re-summation of the leading logarithmic terms in the PS. Actually, some special 
case of the Sudakov form factor can be obtained from T>^(x)as: 



U(Q 2 max ,Q 2 c /(l-x)) 



exp 



V^\x)dx 



(30) 



However, there is one essential difference between T>^ l \x) and the Sudakov form factor in the 
PS: the upper bound of k 2 -^ integration. In the PS, this integration is performed up to the 
energy scale determined by the Born process instead of their maximum value in kinetically 
allowed region. Then, if ull is subtracted in the full phase space of the (A r +l)-particle final 
state, some part of the LL cross section cannot be recovered from the Born cross-section with 
the PS. Then, an appropriate restriction on the phase-space integration must be applied, such 
as 



1 



O'LLsub 



(2 P <l)(2p 2 )v rel Jn vi 



(4) 
N+1 



M 



(4) 
N+1 



hL(x,s)9(Q 2 B -k 2 T ) 



,(31) 



where Q 2 B is the energy scale of the Born process, which must be the same as the factorization 
energy-scale of the PS. 



5 Event generation with a part on shower 
5.1 Conventional method 

The calculations of cross sections for proton- (anti) proton collisions are usually performed as 
follows CHI: 
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Figure 2: Momentum fraction distributions from an x-deterministic parton shower and Cteq5L parton distri- 
bution functions. The distributions of it-quarks (left), u-quarks (center), and gluons (right) are shown. 



1. matrix elements of a given process at a parton level are prepared, 

2. a probability to observe a parton in a proton at given energy scale and momentum fraction 
is obtained from some parton distribution function(PDF) [T^] . 

3. a numerical integration is performed in phase-space with some initial energy distribution 
in the PDF, 

4. for initial-state partons, a (backward) parton shower 20 is applied from the given (high) 
energy scale to a low energy-scale of partons, 

5. for final-state partons, a parton shower is also applied from a given energy scale to close 
to the hadron energy-scale, 

6. some non-perturbative effects according to string fragmentation or color clustering are 
taken into account, and 

7. physical hadrons are formed based on some hadronization model. 

In this method, the parton shower(PS) is used in a supplementary way, just for generating 
multi-partons with a finite transverse momentum. 

5.2 An x-deterministic PS 

The parton shower is a method to solve a DGLAP evolution equation [21] using a Monte Carlo 
method. The PDF's are also obtained by solving the DGLAP equation with experimental 
data-fitting. When an initial distribution of partons is given at some energy scale, the PS 
can reproduce a consistent result as the PDF. For singlet partons, we employ an evolution 
scheme based on momentum distributions rather than those on particle-number distributions, 
which is used for non-singlet processes. The momentum-distribution evolution-scheme has been 
proposed by Tanaka and Munehisa 22 and shows numerically stable results with high efficiency. 

The PDF can give the weight of partons when the momentum fraction (x) and the energy 
scale are given by the users. On the other hand, in the forward evolution scheme of the PS, 
the momentum fraction of a parton (x) is determined only after evolution takes place. This 
method is very inefficient, for instance, for a narrow-resonance production. In order to cure 
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this inefficiency, we have developed an rc-deterministic PS. In a Monte-Carlo procedure in the 
PS, the Q 2 evolution is controlled by the Sudakov form-factor, which gives a non-branching 
probability of partons, and the x determination is done independently according to a splitting 
function. It is not necessary to determine the x value at each branching. After preforming 
all branching procedures using the Sudakov form-factor, the x value at each branching is 
determined to give a total x, being a given value from outside of the PS. In this PS, each event 
has a different weight according to the splitting functions and initial distribution of the PDF. 
When the number of branchings is n, the weight of this event (W) can be obtained as 

1 n 

W = T7rT\P(xi)FpDF(xo)dx l dx , (32) 



n „! 

W o = II / P(xj)dxj, (33) 
j=i 



where P(xi) is a splitting function, Fpdf(xo) is parton distribution function, and xq is a 
momentum fraction before the parton-shower evolution. The momentum fraction after the 
evolution is x = Y\a=o x i- I R order to take this momentum fraction {x) as an independent 
variable given by users, unity 

„ n 

1 = / S(x - x n xj)dx (34) 

J 3=1 

is multiplied to Eq.(|32|) and integration by dx^ is done first as 

1 n / n \ 

— Y[ P(xi)Fp DF (x )6{ x - x Y[ Xj \dxdxidx , (35) 



W 

w °h\ " v- -ji-v 

io = n"^ r - (37) 

llj = l x j 

If we employ an appropriate transformation of the input random numbers, we can expect a 
sufficient efficiency. 

A numerical test of the x-deterministic PS is done by comparing results with the PDF of 
Cteq5L 23 . The initial distributions of partons are taken from Cteq5L at the energy scale of 
the 6-quark mass. Then, after evolution to an energy scale of 100 GeV, the x distributions from 
the PS are compared to those in Cteq5L. One can see in FigEJthe PS reproduced x distributions 
in Cteq5L for both of singlet and non-singlet partons. 

Our procedure to generate QCD events involves the following three steps, : 

2'. numerical integration is performed in the phase-space including x integration, 

3'. the probability for the existence of each parton in the proton at some low-energy scale is 
obtained from the PDF, and 

4'. the ^-deterministic PS performs evolution from a low energy scale to the energy scale of 
a hard scattering while generating multi-partons with a finite transverse momentum, 

instead of steps 2, 3, and 4 in the conventional procedure. 
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Figure 3: Transverse momentum distribution of gluons. The distributions from the PS applied to the Born 
process are shown by * commonly in three histograms. Those distributions from a exac t (left), o exact with double- 
count rejection (middle), and a exac t with double-count rejection and the fey restriction are compared with the 
PS. 



5.3 Negative weight treatment 

Though the LL-subtraction method makes the negative-positive cancellation in the collinear 
region very mild, the NLO amplitude may still have negative values in some region in the 
phase-space. This negative weight is unavoidable due to the perturbative calculation truncated 
at a finite order of expansion. Especially this problem is very serious in QCD due to its large 
expansion coefficient. Even if there are events with negative weight, experimentally measurable 
distributions can be positive with a realistic detector resolution (or a bin-width of histograms). 
Then we decided to employ a 'quasi-unweighted' event generation as follows. 

The BASES 24r| system can treat integrand with negative values and prepare a table of a 
probability density according to the absolute value of the integrand. The SPRING 24| EH] will 
generate events based on that table with unit-weight with sign +1 or — 1 according to a sign 
of integrand at given phase point. This sign will be kept during detector simulations and user 
analysis programs. Finally some bins of histograms will be incremented or decremented by unit 
amount due to their sign. 

6 A test of the LL-subtraction method 

The LL-subtraction method is numerically tested for a process of uu — > /U + /i~ in proton-anti 
proton collision at the CM energy of 2 TeV. Here only non-singlet u-quark is used in this 
test to avoid additional complex problem. Matrix elements of the Born process and the real 
radiation process (uu — > fi + pT gluon) are generated by GRACE. Numerical integration is done 
using a BASES system. A non-singlet u-quark distribution at the energy scale of 4.6 GeV is 
taken from CTEQ5L. The parton evolution from this energy scale to that of hard-scattering, i.e. 
an invariant mass of the muon-pair (s«„), is done using the x-deterministic PS for both of the 
Born and radiative processes. An additional cut of y/s^ > 40 GeV is also applied. 

The distributions of transverse momenta of gluons are shown in Figure 2. For the Born 
process, the largest kx from the PS is filled to histograms if it has kx greater than 1 GeV. For 
the radiative process, that from matrix elements is filled when it has passed the same cut as 
mentioned above. If one does not take care of double-counting of the LL-terms in the matrix 
elements and the PS, the distributions of da exac t/dkx show larger values compared with the 
PS in the Born process, as shown in the left histogram of Figure 2. In order to avoid double- 
counting, we required that kx of the matrix elements be greater than those of any gluons from 
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Figure 4: Energy and transverse momentum distributions of gluons. The * shows the distributions from 
the PS on the Born process, solid histograms from a eX act with double-count rejection, and circles from aLLsub 
combined with those from the PS on the Born process. 

the PS in the radiative process. This double-counting rejection leads to a good agreement of 
the distribution around the low-ftT region, as shown in the middle histogram. However, 
the PS still emits an insufficient number of gluons around a very high-fcr region. If the upper 
bound of the k\ integration is set to be in the calculation with the exact matrix elements, 
fey distributions shows a good agreement, as shown in right histogram. This behavior tells us 
that the parton-shower can reproduce gluon distributions only for the energy scale less than 
s w . This result confirmed the necessity of the step function [Q(Q 2 B — kx)) in Eq.(|31|>. 

The distributions from aiLsub combined with those of the PS on the Born process are 
compared with those from exact matrix-elements with double-counting rejection. In both of 
the gluon energy and the transverse- momentum distributions, the LL-subtr action method can 
give consistent results with those of the exact matrix elements, as shown in Figure 3. 

7 Example:Drell-Yan processes 

In order to demonstrate how our proposed method can really work for physical processes in 
hadron colliders, an example of an event generator for the Drell-Yan process in a proton anti- 
proton (pp) collider of TEVATRON at the next-to-leading order QCD is constructed. All 
possible Feynman diagram contributing to the process pp — > /i + ^ _ (+additional parton) as 
shown in figure 5 are taken into account in the generator. All quarks and gluons except top- 
quark in an initial (anti-)proton are summed up using the PDF of Cteq5L at 6-quark mass of 
4.5 GeV and Q 2 evolution to the energy scale of a hard scattering has been done using the 
x-deterministic parton shower at leading logarithmic order. 

FORTRAN programs to calculate matrix-elements are automatically generated using GRACE 
including loop- and real-radiation diagram as well as tree ones. All quark masses are set to 
be massless in the matrix-element calculations. A QCD coupling constant is calculated us- 
ing a leading-order formula with Aqcd = 146 MeV which is obtained by parameter fittings 
in Cteq5L. The factorization scale is set to be equal to the energy of a ^-pair system. As it 
is known well for the Drell-Yan process, the renormalization energy scale does not explicitly 
appear in the correction terms. Main correction terms are absorbed in the running coupling 
constant and Q 2 evolution of the parton distribution. Remaining terms with the large logarithm 
(L = ln(Q 2 /n)Z) is canceled out after summing up virtual and collinear correction terms. This 
renormalization scale independence is numerically checked in our program. Infrared divergent 
part proportional to \je 2 R and I/em are left in the program and numerically confirmed that 
these divergent coefficients are canceled out after summing up virtual and collinear parts of 
the matrix elements at the order of 10 -20 . No experimental cut except the energy of a /i-pair 
system to be greater than 10 GeV is applied. The total cross section of the Drell-Yan process 
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Figure 5: Feynman diagram of the Drell-Yan process, qq — > /i n + and its radiative processes. 



at 2 TeV of the pp CM energy with above cut is obtained to be 1.026 x 10 3 pb at the tree 
level and 1.288 x 10 3 pb at the NLO, which gives a A'-factor of 1.256. In the cross section 
of 1.288 x 10 3 pb at the NLO, 1.807 x 10 2 pb comes from the visible cross section after the 
LL-subtr action. It is rather good behavior as perturbative correction of QCD. Numerical re- 
sults of the virtual and collinear correction of qq initial state are compared with those based 
on formulae given by Altarelli et al. |26| and are confirmed consistent. It is confirmed that the 
cross sections are independent from the value of Q\ defined in section 4.2 and agree with those 
with simple phase-space slicing method as expected. 

Next let us look at transverse energy distributions of the largest Et jet and a //-pair system. 
Origin of the transverse energy has two sources, one is finite transverse momentum of the PS 
and the other is that from real-radiation matrix-elements. In the LL-subtraction method, 
these two kinds of distributions can be smoothly combined without any double-counting. As 
shown in figure 6, contribution from the parton shower applied on the NLO matrix elements 
(solid histograms) is dominated around small Et region. On the other hand a high Et events 
are dominated by the real-radiation correction (dashed histograms). These two regions are 
smoothly connected as shown by stars (*) which show the sum of solid and dashed histograms. 
If one applied PDF simply on the the real-radiation matrix-elements without any care of double 
counting, small Et events are clearly over estimated as shown by dot-dashed histograms. On 
the other hand, high Et events are consistent with those from LL-subtraction method, because 
they are free from both double-counting and smearing due to the PS. 

8 Conclusions 

A new method to construct event-generators based on next-to-leading order QCD matrix- 
elements and an x-deterministic parton shower is proposed. Matrix elements of loop diagram 
as well as those of a tree level can be generated by the GRACE system. A soft /collinear singularity 
is treated using the leading-log subtraction method. It has been demonstrated that the LL- 
subtraction method can give good agreement with the exact matrix elements without any 
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Figure 6: Transverse energy distributions (-Et) of jets (left figure) and muon-pair system (right figure). For a 
jet distribution, the largest Et among jets from parton showers or matrix elements are filled in the histogram. 
Solid histograms show distributions from parton showers applied on the tree+virtual+collinear collections, 
dashed ones come from real radiation matrix-elements with PS+(LL-subtraction), and dot-dashed ones come 
from the matrix-elements with PDF without any care of the double-counting. 



double-counting problem. The NLO event generator for the Drell-Yan processes has been 
constructed based on our method and showed smooth distributions of transverse momenta 
combining the PS and exact matrix-elements of real radiation processes without any double- 
counting problem. 
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Appendix 



A Effective vertices 
A.l SU(N C ) Color factor 



fabc 

rj-i(2 rpb 



^ y rj-iarj-ia 

a 

2>(T°T 6 ) 

^ ] fabcfbcd 
b,c 



the fundamental representation, 
the structure constants, 

ifabc^- j 

7V 2 -1 



2iV c 



J 



a6 



C G 5 ad = N C S' 



ad 



A. 2 tree vertex 

• quark-gluon vertex: gT a 

• three-gluon vertex: —igfabc 

• four-gluon vertex: —g 2 {f a befcde + cyclic permutation) 



A. 3 quark self-energy 



P 




= nv 2 ) 



on- / off-shell 


self-energy 


p 2 / 
p 2 = 


S(P 2 ) = C F ^(1-In^) 
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A. 4 gluon self-energy 



P 




= n G (p 2 ) 



on- / off-shell 


self-energy 


p 2 = 


neb 2 ) = -c F %(% - § in =£) + r^gcf - | In ^) 
n G (o) = -c F ^^ + r fi iv / ^|_i_ 



A. 5 quark-gluon vertex 




= AJ(pi,p2,fc) 



1 „ v a.s 



A„(pi,P2,fc) = 5T a (C F --C G )^(^V + ^ 



A. 5.1 case I 

k 2 = 0, p 2 = 0,p 2 = q 2 ^ 0, and L = ln=£, where (», j) = (1,2) or (i, j) = (2,1). 
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+ 4L - 10 
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A. 5. 2 case II 

k 2 = q 2 , v\=v\ = 0, and L = In 
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2 2L-A Q 
€IR 


+ 4 + 3L-L 2 


H 







H 1 


£IR 




•r 2 








A. 6 three gluon vertex 



0.3 




k 2 = k% = 0,k1 = q 2 ^0,L = ln- 



I 1 



Qi 



-iqf aia2a3 — 
4vr 



3=1 



£/i:gluon loop, (^ghost loop, C?3:gluon loop (fish type), ^iquark loop 
Independent momenta 

• v{ xmWi = (h - k 2 Y 3 g^ 2 



P. 
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Cij 


3 = 1 


J = 2 


J = 3 


i = 1 


4~ R 5+L) 

-^(103 - 5LL + 9L 2 ) + ^ 
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